Linear Wave 방정식#
강좌: 기초 전산유체역학
선형 파동 방정식#
1차원 선형 파동 방정식은 다음과 같다.
여기서 상수 \(a\) 는 파의 전파속도이다. 다음 초기 조건에 의한 완전해는
아래와 같다.
유한차분법#
\(n\) 번째 시간 \(t^n\), \(j\) 번째 격자점 \(x_j\) 에서 근사해를 다음과 같이 표현하자.
First approach (Central difference)#
시간에 대한 차분은 우선 Euler Explicit 방법을 생각하자.
공간에 대한 차분은 중앙차분을 생각하자.
이를 정리하면 다음 기법을 만들 수 있다.
Second approach (Upwind difference)#
\(a>0\) 인 경우 공간에 대한 차분을 1차 backward difference로 생각하자.
이를 정리하면 다음 기법을 만들 수 있다.
Sine wave 예제#
계산 영역은 \(x \in [-1, 1]\) 이고 초기 조건은 다음과 같다.
양 끝점에서 경계 조건은 Periodic 조건을 주자. 수식적으로는
시간 \(t=1.5\) 일 때 해를 구하자.
계산 격자 구성, Solution array 구성#
계산 영역을 \(n_x + 1\)개의 점으로 나누어보자.
즉 격자점은 다음과 같다.

Fig. 8 Grid#
첫번째 격자점을 계산하기 위해서는 \(x_0 = -1 - \Delta x\) 값이 필요하다. 또한 마지막 격자점을 계산하기 위해서는 \(x_{n_x+2} = 1 + \Delta x\) 값이 필요하다.
이를 대칭조건으로 구현하면 다음과 같다.
경계 조건을 위해서 Solution array는 격자점 + 2개의 경계 조건을 고려해서 \(n_x+3\) 으로 구성한다.
시간 간격은 \(\Delta t = 0.01\) 에 대해서 \(n_x=50\) 인 경우 계산하면 다음과 같다.
from matplotlib import pyplot as plt
import numpy as np
plt.style.use('ggplot')
plt.rcParams['figure.dpi'] = 150
def central(nx, u, dt, dx, a, du):
"""
Central difference
Parameters
----------
nx : integer
등분 개수
u : array
Solution
dt : float
시간간격
dx : float
격자 간격
a : float
Wave speed
du : array
증가량
"""
for i in range(1, nx+2):
du[i] = -0.5*a*(u[i+1] - u[i-1])/dx*dt
def central_v1(u, dt, dx, a):
"""
Central difference (vector version)
Parameters
----------
u : array
Solution
dt : float
시간간격
dx : float
격자 간격
a : float
Wave speed
Return
------
du : array
증가량
"""
# Size of u : nx + 3
# 0, nx+2 (or -1) - ghost cell,
# 1 - nx+1 (or -2) - Inner
return -0.5*a*(u[2:] - u[:-2])/dx*dt
def bc_periodic(u):
"""
Boundary condition (peridoic)
Parameter
---------
u : array
solution
"""
# index (nx : -3), (nx+2 : -1)
u[0] = u[-3]
u[-1] = u[2]
# Constants
a = 1.0
nx = 50
dt = 0.01
t_target = 1.5
# Make grid
x = np.linspace(-1, 1, nx+1)
dx = np.diff(x)[0]
# Solution array
u = np.empty(nx+3)
du = np.zeros_like(u)
# Initialize
u[1:-1] = np.sin(np.pi*x)
bc_periodic(u)
# Calculation
t = 0
while abs(t - t_target) > 1e-8:
# Adjust time step to reach target time
dt = min(dt, t_target - t)
# Periodic BC
bc_periodic(u)
# Scheme
#du[1:-1] = central_v1(u, dt, dx, a)
central(nx, u, dt, dx, a, du)
# Update
u += du
t += dt
# Compare with exact solution
u_exact = np.sin(np.pi*(x-a*t))
plt.plot(x, u_exact)
plt.plot(x, u[1:-1])
plt.legend(['Exact', 'Computed'])
<matplotlib.legend.Legend at 0xff3314c95d00>

실습#
Central 기법에 대해 격자 개수와 시간 간격을 아래와 같이 달리하면셔 결과를 비교하시오.
\(N\) = 50, 100, 200
\(\Delta t = 0.005, 0.01, 0.02, 0.05, 0.1\)
Upwind difference를 구현한 후, Sine Wave 문제를 해석하시오. Central 기법과 결과 차이를 비교하시오.
다음 Square Wave에 대해 Central difference 와 Upwind difference를 해석하시오.
(Optional) \(a < 0\) 일때 Upwind difference 기법은 공간에 대해 Forward 차분식을 적용해야 한다. \(a=-1\)인 경우에 대해 Upwind 기법을 구현하시오.